package com.informatics.polymer.server.Polymerisation;

import java.io.IOException;
import java.util.ArrayList;
import java.util.List;
import java.util.Random;

import nu.xom.Builder;
import nu.xom.Document;
import nu.xom.ParsingException;
import nu.xom.ValidityException;

import org.xmlcml.cml.element.CMLAtom;
import org.xmlcml.cml.element.CMLBond;
import org.xmlcml.cml.element.CMLFragment;
import org.xmlcml.cml.element.CMLJoin;
import org.xmlcml.cml.element.CMLMolecule;

/**
 * Class to create random copolymers.<br>
 * </br> <b>Example</b><br>
 * </br> RandomCopolymerisation randomCopolymerisation = new
 * RandomCopolymerisation
 * (0.9,0.2,0.01,0.05,0.8,0.01,180.0,"[R]c1cccc(C[R])c1","[R]OCC1CCCC(C[R])C1"
 * ,"CC[R]","CCO[R]");<br>
 * </br> String randomCopolymer2D =randomCopolymerisation.getPolymer().toXML();<br>
 * </br> String randomCopolymer3D =
 * randomCopolymerisation.add3DCoordinates(randomCopolymer2D);<br>
 * </br> <b>Note</b><br>
 * </br> Coordinate generation is done using external dependencies, from the
 * Chemistry Development kit<br>
 * </br> <a href="http://pele.farmbio.uu.se/nightly/api/org/openscience/cdk/modeling/builder3d/ModelBuilder3D.html"
 * >ModelBuilder3D class</a>.
 * 
 * @author Ed Cannon
 * @version 1.0
 */
public class RandomCopolymerisation implements Polymerisation {
	/**
	 * Maximum torsion angle allowed.
	 */
	private static final int TORSION_ANGLE_UPPER_LIMIT = 360;
	/**
	 * Double value representing the probability of repeat unit A occurring
	 * after A.
	 */
	private double pAA;
	/**
	 * Double value representing the probability of repeat unit B occurring
	 * after A.
	 */
	private double pAB;
	/**
	 * Double value representing the probability of termination after adding A.
	 */
	private double pTA;
	/**
	 * Double value representing the probability of repeat unit A occurring
	 * after B.
	 */
	private double pBA;
	/**
	 * Double value representing the probability of repeat unit B occurring
	 * after B.
	 */
	private double pBB;
	/**
	 * Double value representing the probability of termination after adding B.
	 */
	private double pTB;
	/**
	 * Double value representing the angle between repeat units.
	 */
	private double torsionAngle, x;
	/**
	 * Double value representing the rate of addition of repeat unit A.
	 */
	private double dA;
	/**
	 * Double value representing the rate of addition of repeat unit B.
	 */
	private double dB;
	/**
	 * Integer value representing the maximum number of repeat unit As in the
	 * polymer.
	 */
	private int nA;
	/**
	 * Integer value representing the maximum number of repeat unit Bs in the
	 * polymer.
	 */
	private int nB;
	/**
	 * Boolean to test if input probability is correct.
	 */
	private boolean val;
	/**
	 * CMLMolecule representing the randomCopolymer.
	 */
	private CMLMolecule randomCopolymer;
	/**
	 * String reprsenting the CMLMolecule to add.
	 */
	private String fragmentValue;
	/**
	 * CMLMolecule representation of repeat unit A.
	 */
	private CMLMolecule repeatUnitA;
	/**
	 * CMLMolecule representation of repeat unit B.
	 */
	private CMLMolecule repeatUnitB;
	/**
	 * CMLMolecule representation of end group 1.
	 */
	private CMLMolecule endGroup1;
	/**
	 * CMLMolecule representation of end group 2.
	 */
	private CMLMolecule endGroup2;

	/**
	 * ArrayList containing all ids for building the polymer.
	 */
	private ArrayList<String> ids;
	/**
	 * Constructs an object to undertake random copolymerisation.
	 */
	public RandomCopolymerisation() {
	}

	/**
	 * Constructs an object to undertake random copolymerisation using the
	 * probability of repeat unit A occurring after A, probability of a repeat
	 * unit B coming after a repeat unit A, the probability of terminating the
	 * polymerisation after adding repeat unit A, the probability of a repeat
	 * unit A coming after a repeat unit B, the probability of repeat unit B
	 * occurring after B, the probability of terminating the polymerisation
	 * after repeat unit B, the torsion angle and SMILE strings for end groups 1
	 * and 2 and repeat units A and B.
	 * 
	 * @param pAAValue
	 *            - double value representing the probability of a repeat unit A
	 *            coming after a repeat unit A.
	 * @param pABValue
	 *            - double value representing the probability of a repeat unit B
	 *            coming after a repeat unit A.
	 * @param pTAValue
	 *            - double value representing the probability of terminating the
	 *            polymerisation after a repeat unit A is added to the chain.
	 * @param pBAValue
	 *            - double value representing the probability of a repeat unit A
	 *            coming after a repeat unit B.
	 * @param pBBValue
	 *            - double value representing the probability of a repeat unit B
	 *            coming after a repeat unit B.
	 * @param pTBValue
	 *            - double value representing the probability of terminating the
	 *            polymerisation after a repeat unit B is added to the chain.
	 * @param torsionAngleValue
	 *            - double value representing the torsion angle between
	 *            subsequent repeat units.
	 * @param ruAValue
	 *            - SMILE string representation of repeat unit A.
	 * @param ruBValue
	 *            - SMILE string representation of repeat unit B.
	 * @param eg1Value
	 *            - SMILE string representation of end group 1.
	 * @param eg2Value
	 *            - SMILE string representation of end group 2.
	 */
	public RandomCopolymerisation(final double pAAValue, final double pABValue,
			final double pTAValue, final double pBAValue,
			final double pBBValue, final double pTBValue,
			final double torsionAngleValue, final String ruAValue,
			final String ruBValue, final String eg1Value, final String eg2Value) {
		checkParameterSetup(pAAValue, pABValue, pTAValue, pBAValue, pBBValue,
				pTBValue, torsionAngleValue);
		checkStringMoleculeSetup(ruAValue, ruBValue, eg1Value, eg2Value);
		setUpComponents();
	}

	/**
	 * Constructs an object to undertake random copolymerisation using the
	 * probability of repeat unit A occurring after A, probability of a repeat
	 * unit B coming after a repeat unit A, the probability of terminating the
	 * polymerisation after adding repeat unit A, the probability of a repeat
	 * unit A coming after a repeat unit B, the probability of repeat unit B
	 * occurring after B, the probability of terminating the polymerisation
	 * after repeat unit B, the torsion angle and CMLMolecules for end groups 1
	 * and 2 and repeat units A and B.
	 * 
	 * @param pAAValue
	 *            - double value representing the probability of a repeat unit A
	 *            coming after a repeat unit A.
	 * @param pABValue
	 *            - double value representing the probability of a repeat unit B
	 *            coming after a repeat unit A.
	 * @param pTAValue
	 *            - double value representing the probability of terminating the
	 *            polymerisation after a repeat unit A is added to the chain.
	 * @param pBAValue
	 *            - double value representing the probability of a repeat unit A
	 *            coming after a repeat unit B.
	 * @param pBBValue
	 *            - double value representing the probability of a repeat unit B
	 *            coming after a repeat unit B.
	 * @param pTBValue
	 *            - double value representing the probability of terminating the
	 *            polymerisation after a repeat unit B is added to the chain.
	 * @param torsionAngleValue
	 *            - double value representing the torsion angle between
	 *            subsequent repeat units.
	 * @param repeatUnitAValue
	 *            - CMLMolecule representation of repeat unit A.
	 * @param repeatUnitBValue
	 *            - CMLMolecule representation of repeat unit B.
	 * @param endGroup1Value
	 *            - CMLMolecule representation of end group 1.
	 * @param endGroup2Value
	 *            - CMLMolecule representation of end group 2.
	 */
	public RandomCopolymerisation(final double pAAValue, final double pABValue,
			final double pTAValue, final double pBAValue,
			final double pBBValue, final double pTBValue,
			final double torsionAngleValue, final CMLMolecule repeatUnitAValue,
			final CMLMolecule repeatUnitBValue,
			final CMLMolecule endGroup1Value, final CMLMolecule endGroup2Value) {
		checkParameterSetup(pAAValue, pABValue, pTAValue, pBAValue, pBBValue,
				pTBValue, torsionAngleValue);
		checkCMLMoleculeSetUp(repeatUnitAValue, repeatUnitBValue,
				endGroup1Value, endGroup2Value);
		setUpComponents();
	}

	/**
	 * Sets the probability of repeat unit A occurring after repeat unit A.
	 * 
	 * @param pAAValue
	 *            - double value representing the probability of a repeat unit A
	 *            coming after a repeat unit A
	 */
	public final void setProbabilityOfAAfterA(final double pAAValue) {
		this.pAA = pAAValue;
	}

	/**
	 * Sets the probability of repeat unit B occurring after repeat unit A.
	 * 
	 * @param pABValue
	 *            - double value representing the probability of a repeat unit A
	 *            coming after a repeat unit B
	 */
	public final void setProbabilityOfBAfterA(final double pABValue) {
		this.pAB = pABValue;
	}

	/**
	 * Sets the probability of termination after repeat unit A is added.
	 * 
	 * @param pTAValue
	 *            - double value representing the probability of terminating the
	 *            polymerisation after a repeat unit A is added to the chain
	 */
	public final void setProbabilityOfATermination(final double pTAValue) {
		this.pTA = pTAValue;
	}

	/**
	 * Sets the probability of repeat unit A occurring after repeat unit B.
	 * 
	 * @param pBAValue
	 *            - double value representing the probability of a repeat unit B
	 *            coming after a repeat unit A
	 */
	public final void setProbabilityOfAAfterB(final double pBAValue) {
		this.pBA = pBAValue;
	}

	/**
	 * Sets the probability of termination after repeat unit B is added.
	 * 
	 * @param pTBValue
	 *            - double value representing the probability of terminating the
	 *            polymerisation after a repeat unit B is added to the chain
	 */
	public final void setProbabilityOfBTermination(final double pTBValue) {
		this.pTB = pTBValue;
	}

	/**
	 * Sets the probability of repeat unit B occurring after repeat unit B.
	 * 
	 * @param pBBValue
	 *            double value representing the probability of a repeat unit B
	 *            coming after a repeat unit B
	 */
	public final void setProbabilityOfBAfterB(final double pBBValue) {
		this.pBB = pBBValue;
	}

	/**
	 * Sets the torsion angle between repeat units.
	 * 
	 * @param torsionAngleValue
	 *            double value representing the torsion angle between subsequent
	 *            repeat units
	 */
	public final void setTorsionAngle(final double torsionAngleValue) {
		this.torsionAngle = torsionAngleValue;
	}

	/**
	 * Sets CMLMolecule to be repeat unit A.
	 * 
	 * @param repeatUnitAMolecule
	 *            - CMLMolecule representation of repeat unit A
	 */
	public final void setRepeatUnitA(final CMLMolecule repeatUnitAMolecule) {
		this.repeatUnitA = repeatUnitAMolecule;
	}

	/**
	 * Sets CMLMolecule to be repeat unit B.
	 * 
	 * @param repeatUnitBMolecule
	 *            - CMLMolecule representation of repeat unit B
	 */
	public final void setRepeatUnitB(final CMLMolecule repeatUnitBMolecule) {
		this.repeatUnitB = repeatUnitBMolecule;
	}

	/**
	 * Sets CMLMolecule to be end group 1.
	 * 
	 * @param eg1
	 *            - CMLMolecule representation of end group 1
	 */
	public final void setEndGroup1(final CMLMolecule eg1) {
		this.endGroup1 = eg1;
	}

	/**
	 * Sets CMLMolecule to be end group 2.
	 * 
	 * @param eg2
	 *            - CMLMolecule representation of end group 2
	 */
	public final void setEndGroup2(final CMLMolecule eg2) {
		this.endGroup2 = eg2;
	}

	/**
	 * Returns a double value representing the probability of a repeat unit A
	 * coming after a repeat unit A.
	 * 
	 * @return pAA
	 */
	// Have option to get all values
	public final double getProbabilityOfAA() {
		return pAA;
	}

	/**
	 * Returns a double value representing the probability of a repeat unit B
	 * coming after a repeat unit A.
	 * 
	 * @return pAB
	 */
	public final double getProbabilityOfAB() {
		return pAB;
	}

	/**
	 * Returns a double value representing the probability of terminating the
	 * polymerisation after a repeat unit A is added to the chain.
	 * 
	 * @return pTA
	 */
	public final double getProbabilityOfPTA() {
		return pTA;
	}

	/**
	 * Returns a double value representing the probability of a repeat unit A
	 * coming after a repeat unit B.
	 * 
	 * @return pBA
	 */
	public final double getProbabilityofBA() {
		return pBA;
	}

	/**
	 * Returns a double value representing the probability of a repeat unit B
	 * coming after a repeat unit B.
	 * 
	 * @return pBB
	 */
	public final double getProbabilityofBB() {
		return pBB;
	}

	/**
	 * Returns a double value representing the probability of terminating the
	 * polymerisation after a repeat unit B is added to the chain.
	 * 
	 * @return pTB
	 */
	public final double getProbabilityofPTB() {
		return pTB;
	}

	/**
	 * Returns a double value representing the torsion angle between subsequent
	 * repeat units.
	 * 
	 * @return torsionAngle
	 */
	public final double getTorsionAngle() {
		return torsionAngle;
	}

	/**
	 * Returns a CMLMolecule representation of repeat unit A.
	 * 
	 * @return repeatUnitA
	 */
	public final CMLMolecule getRepeatUnitA() {
		return repeatUnitA;
	}

	/**
	 * Returns a CMLMolecule representation of repeat unit B.
	 * 
	 * @return repeatUnitB
	 */
	public final CMLMolecule getRepeatUnitB() {
		return repeatUnitB;
	}

	/**
	 * Returns a CMLMolecule representation of end group 1.
	 * 
	 * @return endGroup1
	 */
	public final CMLMolecule getEndGroup1() {
		return endGroup1;
	}

	/**
	 * Returns a CMLMolecule representation of end group 2.
	 * 
	 * @return endGroup2
	 */
	public final CMLMolecule getEndGroup2() {
		return endGroup2;
	}

	/**
	 * Returns false if probability is within the range of 0 - 1 else true if
	 * the probability is not correct.
	 * 
	 * @param probabilityValue
	 *            - probability value to test if in the range of 0 - 1
	 * @throws tag
	 *             for 'IllegalArgumentException
	 * @return boolean probability check.
	 */
	private boolean testProbability(final double probabilityValue)
			throws IllegalArgumentException {
		boolean failed = false;
		if (probabilityValue < 0 || probabilityValue > 1) {
			failed = true;
			throw new IllegalArgumentException(
					"A probability must be in the range of 0-1");
		} // else {
		// failed = false;
		// }
		return failed;
	}

	/**
	 * Sets the maximum number of repeat unit A components in the polymer.
	 */
	private void setLengthOfA() {
		nA = (int) (1 / pBA);
	}

	/**
	 * Sets the maximum number of repeat unit B components in the polymer.
	 */
	private void setLengthOfB() {
		nB = (int) (1 / pAB);
	}

	/**
	 * Updates the mol fractions when a new repeat unit is added.
	 */
	private void updateMolFractions() {
		try {
			x = nA / nB;
		} catch (ArithmeticException e) {
			e.printStackTrace();
		}
	}

	/**
	 * Checks which repeat unit to add.
	 */
	private void updateRateValues() {
		if (val == true) { // then add after A, else add after b probability
			dA = 1 + (pAA * x);
			dB = 1 + (pBA / x);
		} else {
			dA = 1 + (pAB * x);
			dB = 1 + (pBB / x);
		}
	}

	/*
	 * Returns fragment value.
	 * 
	 * @return fragmentValue
	 */
	/*
	 * private String getFragmentValue() { return fragmentValue; }
	 */
	/**
	 * Sets next fragment in random copolymer.
	 */
	private void nextFragment() {
		fragmentValue = null;
		updateRateValues();
		double result = dA / dB;
		if (result > 1) {
			double aOrTerminate = dA / pTA;
			if (aOrTerminate > 1) {
				fragmentValue = "RuA";
				nA--;
			} else {
				fragmentValue = "Eg2";
			}
		} else if (result < 1) {
			double bOrTerminate = dB / pTB;
			if (bOrTerminate > 1) {
				fragmentValue = "RuB";
				nB--;
			} else {
				fragmentValue = "Eg2";
			}
		} else if (result == 1) {
			// Randomly choose a new number
			Random randomGenerator = new Random();
			int randomInt = randomGenerator.nextInt(2);
			if (randomInt == 0) {
				fragmentValue = "RuA";
				nA--;
			} else {
				fragmentValue = "RuB";
				nB--;
			}
		}
	}

	/**
	 * Returns a list of ids of repeat units and their order in which they
	 * should be added in order to build the polymer.
	 * 
	 * @return ArrayList<String>
	 */
	private ArrayList<String> randomCopolymerSeparateComponents() {
		setLengthOfA();
		setLengthOfB();
		updateMolFractions(); // just use molFractionA in calculation to be
		// consistent with Google Book
		ArrayList<String> repeatUnitsToAdd = new ArrayList<String>();
		if (fragmentValue == null) {
			nextFragment();
		}
		repeatUnitsToAdd.add(fragmentValue);
		int totalCounter = nA + nB;
		int counter = 0;
		while (fragmentValue != "Eg2") {
			updateMolFractions();
			nextFragment();
			repeatUnitsToAdd.add(fragmentValue);
			counter++;
			if (counter == totalCounter - 1) {
				break;
			}
		}
		repeatUnitsToAdd.add("Eg2");
		repeatUnitsToAdd.add(0, "Eg1");
		return repeatUnitsToAdd;
	}

	// public void setVal(boolean b){
	// val = b;
	// }
	/**
	 * Returns an ArrayList<CMLMolecule> based on the pre-generated repeat unit
	 * string ids.
	 * 
	 * @param ids
	 *            - ArrayList<String> representing the fragments to be added in
	 *            order to build the CMLMolecule fragments for the polymer
	 * @return ArrayList<CMLMolecule>
	 */
	private ArrayList<CMLMolecule> getMolecules(final ArrayList<String> ids) {
		int repeatUnitACounter = 0;
		int repeatUnitBCounter = 0;
		ArrayList<CMLMolecule> moleculeList = new ArrayList<CMLMolecule>();
		moleculeList.add(endGroup1);
		// Construct ArrayList
		for (int i = 0; i < ids.size(); i++) {
			if (ids.get(i).equals("RuA")) {
				CMLMolecule ruAMolecule = countExpression(repeatUnitACounter,
						repeatUnitA, "RuA");
				moleculeList.add(ruAMolecule);
				repeatUnitACounter++;
			} else if (ids.get(i).equals("RuB")) {
				CMLMolecule ruBMolecule = countExpression(repeatUnitBCounter,
						repeatUnitB, "RuB");
				moleculeList.add(ruBMolecule);
				repeatUnitBCounter++;
			} else if (ids.get(i).equals("Eg2")) {
				moleculeList.add(endGroup2);
			}
		}
		return moleculeList;
	}

	/**
	 * Polymerises repeat unit according to counter Value.
	 * 
	 * @param counterValue
	 *            - int parameter representing the number of repeat units to add
	 *            of a specific repeat unit type
	 * @param ru
	 *            - CMLMolecule representation of the repeat unit
	 * @param extraId
	 *            - String representation of a unique id assigned to the repeat
	 *            unit
	 * @return CMLMolecule
	 */
	private CMLMolecule countExpression(final int counterValue,
			final CMLMolecule ru, final String extraId) {
		CMLMolecule repeatUnitWithCountExpression = new CMLMolecule();
		for (int i = 0; i < ru.getAtomCount(); i++) {
			CMLAtom a = new CMLAtom();
			a.setAttribute("id", extraId + Integer.toString(counterValue)
					+ ru.getAtom(i).getId());
			a.setAttribute("elementType", ru.getAtom(i).getElementType());
			a.setFormalCharge(ru.getAtom(i).getFormalCharge());
			// Add aromaticity flag
			if (ru.getAtom(i).getAttribute("aromatic") != null) {
				a.setAttribute("aromatic", ru.getAtom(i).getAttributeValue(
						"aromatic"));
			}
			repeatUnitWithCountExpression.addAtom(a);
		}
		for (int i = 0; i < ru.getBondCount(); i++) {
			CMLBond b = new CMLBond();
			String id1 = ru.getBonds().get(i).getAtomId(0);
			String id2 = ru.getBonds().get(i).getAtomId(1);
			b.setAttribute("atomRefs2", extraId
					+ Integer.toString(counterValue) + id1 + " " + extraId
					+ Integer.toString(counterValue) + id2);
			b.setId(extraId + Integer.toString(counterValue) + id1 + " "
					+ extraId + Integer.toString(counterValue) + id2);
			b.setOrder(ru.getBonds().get(i).getOrder());
			repeatUnitWithCountExpression.addBond(b);
		}
		return repeatUnitWithCountExpression;
	}

	/**
	 * Returns and XML string representation of the generated random copolymer
	 * in 2D.
	 * 
	 * @return randomCopolymer.toXML()
	 */
	public final String generatePML() {
		return randomCopolymer.toXML();
	}

	/**
	 * . Method used to clean up molecule by eliminating R-R bonds
	 */
	private void eliminateRRBonds() {
		List<CMLBond> allBonds = randomCopolymer.getBonds();
		ArrayList<CMLBond> bothRAtoms2 = new ArrayList<CMLBond>();
		for (int i = 0; i < allBonds.size(); i++) {
			if ("R".equals(allBonds.get(i).getAtom(0).getElementType())
					&& ("R".equals(allBonds.get(i).getAtom(1).getElementType()))) {
				bothRAtoms2.add(allBonds.get(i));
			}
		}
		for (int r = 0; r < bothRAtoms2.size(); r++) {
			RRBonds rrval = new RRBonds(bothRAtoms2.get(r), randomCopolymer);
			rrval.automaticallySetX1X2();
			rrval.createX1X2Bond();
			rrval.deleteRgroupsAndBonds();
			randomCopolymer = rrval.getUpdatedMolecule();
		}
	}

	/**
	 * Adds 3D coordinates to a pre-built polymer representation in 2D.
	 * 
	 * @param pmlStringNoCoordinates
	 *            - 2D XML representation of the built polymer
	 * @return randomCopolymer with 3D coordinates.
	 * @throws Exception
	 *             ModelBuilder3D class @link
	 *             (http://pele.farmbio.uu.se/nightly/
	 *             api/org/openscience/cdk/modeling
	 *             /builder3d/ModelBuilder3D.html).
	 * @throws Exception
	 *             if 3D coordinate generation fails.
	 */
	public final String add3DCoordinates(final String pmlStringNoCoordinates)
			throws Exception {
		String molWithCoords = PolymerisationUtils
				.add3DCoordinates(pmlStringNoCoordinates);
		return molWithCoords;
	}

	/*
	 * public static void main(String[] args) throws Exception {
	 * RandomCopolymerisation randomCopolymerisation = new
	 * RandomCopolymerisation
	 * (0.9,0.2,0.01,0.05,0.8,0.01,180.0,"[R]c1cccc(C[R])c1"
	 * ,"[R]OCC1CCCC(C[R])C1","CC[R]","CCO[R]"); CMLFragment frag =
	 * randomCopolymerisation.getFragment(); Document pmlBasicDoc
	 * =PolymerisationUtils.createPMLBasic(frag); Document pmlIntermediateDoc
	 * =PolymerisationUtils.createPMLIntermediate(frag);
	 * System.out.println(pmlBasicDoc.toXML());
	 * System.out.println(pmlIntermediateDoc.toXML()); } ///Document pmlBasicDoc
	 * = / randomCopolymerisation.createPMLBasic(); Document pmlIntermediateDoc
	 * = randomCopolymerisation.createPMLIntermediate();
	 * System.out.println(pmlBasicDoc.toXML());
	 * System.out.println(pmlIntermediateDoc.toXML()); // String
	 * mol=randomCopolymerisation.getPolymer().toXML();
	 * //System.out.println(mol); //String molWithcoords =
	 * randomCopolymerisation.add3DCoordinates(mol); //
	 * //System.out.println(molWithcoords); //PrintWriter out= new
	 * PrintWriter(new // BufferedWriter(new FileWriter("C://foo.cml"),16384));
	 * // //out.write(molWithcoords); //Assume start with repeat unit A // } }
	 */
	/**
	 * Returns CMLMolecule representation of the built polymer in 2D.
	 * 
	 * @return randomCopolymer
	 */
	public final CMLMolecule getPolymer() {
		//boolean buildable = PolymerisationUtils.isBuildable(getFragment());
		boolean buildable = true;
		if (buildable == true) {
			randomCopolymer = new CMLMolecule();
			// ArrayList<String> ids = randomCopolymerSeparateComponents();
			ArrayList<CMLMolecule> m = getMolecules(ids);
			ArrayList<ArrayList<String>> rGroups = PolymerisationUtils
					.getRgroupConnectors(m);
			for (int i = 0; i < m.size(); i++) {
				for (int j = 0; j < m.get(i).getAtomCount(); j++) {
					CMLAtom a = new CMLAtom();
					a.setAttribute("id", m.get(i).getAtom(j).getId());
					a.setAttribute("elementType", m.get(i).getAtom(j)
							.getElementType());
					a.setFormalCharge(m.get(i).getAtom(j).getFormalCharge());
					// Add aromaticity flag
					if (m.get(i).getAtom(j).getAttribute("aromatic") != null) {
						a.setAttribute("aromatic", m.get(i).getAtom(j)
								.getAttributeValue("aromatic"));
					}
					randomCopolymer.addAtom(a);
				}
				for (int j = 0; j < m.get(i).getBondCount(); j++) {
					CMLBond b = new CMLBond();
					String id1 = m.get(i).getBonds().get(j).getAtomId(0);
					String id2 = m.get(i).getBonds().get(j).getAtomId(1);
					b.setAttribute("atomRefs2", id1 + " " + id2);
					b.setId(id1 + " " + id2);
					b.setOrder(m.get(i).getBonds().get(j).getOrder());
					randomCopolymer.addBond(b);
				}
			}
			// Add R-R bonds
			for (int i = 0; i < rGroups.size() - 1; i++) {
				if (i == 0) {
					CMLBond b = new CMLBond();
					b.setAttribute("atomRefs2", rGroups.get(0).get(0) + " "
							+ rGroups.get(1).get(0));
					b
							.setId(rGroups.get(0).get(0) + " "
									+ rGroups.get(1).get(0));
					b.setOrder("1");
					randomCopolymer.addBond(b);
				} else if (i > 0 && i < rGroups.size() - 1 && (i % 2 == 0)) {
					CMLBond b = new CMLBond();
					b.setAttribute("atomRefs2", rGroups.get(i).get(1) + " "
							+ rGroups.get(i + 1).get(0));
					b.setId(rGroups.get(i).get(1) + " "
							+ rGroups.get(i + 1).get(0));
					b.setOrder("1");
					randomCopolymer.addBond(b);
				} else if (i > 0 && i < rGroups.size() - 1 && (i % 2 != 0)) {
					CMLBond b = new CMLBond();
					b.setAttribute("atomRefs2", rGroups.get(i).get(1) + " "
							+ rGroups.get(i + 1).get(0));
					b.setId(rGroups.get(i).get(1) + " "
							+ rGroups.get(i + 1).get(0));
					b.setOrder("1");
					randomCopolymer.addBond(b);
				}

				if (i == rGroups.size() - 2) {
					CMLBond b = new CMLBond();
					b.setAttribute("atomRefs2", rGroups.get(rGroups.size() - 2)
							.get(1)
							+ " " + rGroups.get(rGroups.size() - 1).get(0));
				}
			}
			// Remove R-R bonds
			eliminateRRBonds();
			return randomCopolymer;
		} else {
			UnbuildableException unbuildableException = new UnbuildableException(
					"Can't be built");
			unbuildableException.getMessage();
			return null;
		}
	}

	/**
	 * Checks parameter set up for the constructor is fine.
	 * 
	 * @param pAAValue
	 *            - double representing the probability of repeat unit A after
	 *            A.
	 * @param pABValue
	 *            - double representing the probability of repeat unit B after
	 *            A.
	 * @param pTAValue
	 *            - double representing the probability of termination after A.
	 * @param pBAValue
	 *            - double representing the probability of repeat unit A after
	 *            B.
	 * @param pBBValue
	 *            - double representing the probability of repeat unit B after
	 *            B.
	 * @param pTBValue
	 *            - double representing the probability of termination after B.
	 * @param torsionAngleValue
	 *            - double representing the angle between repeat units.
	 */
	private void checkParameterSetup(final double pAAValue,
			final double pABValue, final double pTAValue,
			final double pBAValue, final double pBBValue,
			final double pTBValue, final double torsionAngleValue) {
		boolean valProb = testProbability(pAAValue);
		if (valProb == false) {
			this.pAA = pAAValue;
		}
		boolean valAB = testProbability(pABValue);
		if (valAB == false) {
			this.pAB = pABValue;
		}
		boolean valTA = testProbability(pTAValue);
		if (valTA == false) {
			this.pTA = pTAValue;
		}
		boolean valBA = testProbability(pBAValue);
		if (valBA == false) {
			this.pBA = pBAValue;
		}
		boolean valBB = testProbability(pBBValue);
		if (valBB == false) {
			this.pBB = pBBValue;
		}
		boolean valTB = testProbability(pTBValue);
		if (valTB == false) {
			this.pTB = pTBValue;
		}
		if (torsionAngleValue < 0) {
			throw new IllegalArgumentException("Negative torsion angle!");
		} else if (torsionAngleValue > TORSION_ANGLE_UPPER_LIMIT) {
			this.torsionAngle = Math.round(torsionAngleValue
					/ TORSION_ANGLE_UPPER_LIMIT);
		} else {
			this.torsionAngle = torsionAngleValue;
		}

	}

	/**
	 * Checks to constructor with SMILE string representation is ok.
	 * 
	 * @param ruA
	 *            - SMILE string representation of repeat unit A.
	 * @param ruB
	 *            - SMILE string representation of repeat unit B.
	 * @param eg1
	 *            - SMILE string representation of end group 1.
	 * @param eg2
	 *            - SMILE string representation of end group 2.
	 */
	private void checkStringMoleculeSetup(final String ruA, final String ruB,
			final String eg1, final String eg2) {
		if (ruA == null) {
			throw new NullPointerException(
					"No repeat unit entered, please enter a repeat unit smile string");
		} else {
			int rCount = PolymerisationUtils.countNumberofFragmentPoints(ruA,
					"[R]");
			if (rCount == 2) {
				this.repeatUnitA = PolymerisationUtils.convertSmileString2CML(
						ruA, "RUA");
			} else if (rCount != 2) {
				throw new IllegalArgumentException(
						"Repeat units must have two R-groups!");
			}
		}

		if (ruB == null) {
			throw new NullPointerException(
					"No repeat unit entered, please enter a repeat unit smile string");
		} else {
			int rCount = PolymerisationUtils.countNumberofFragmentPoints(ruB,
					"[R]");
			if (rCount == 2) {
				this.repeatUnitB = PolymerisationUtils.convertSmileString2CML(
						ruB, "RUB");
			} else if (rCount != 2) {
				throw new IllegalArgumentException(
						"Repeat units must have two R-groups!");
			}
		}
		if (eg1 == null
				|| PolymerisationUtils.countNumberofFragmentPoints(eg1, "[R]") != 1) {
			this.endGroup1 = PolymerisationUtils.convertSmileString2CML(
					"[R][H]", "Eg1");
		} else {
			int rCount = PolymerisationUtils.countNumberofFragmentPoints(eg1,
					"[R]");
			if (rCount == 1) {
				this.endGroup1 = PolymerisationUtils.convertSmileString2CML(
						eg1, "Eg1");
			}
		}
		if (eg2 == null) {
			this.endGroup2 = PolymerisationUtils.convertSmileString2CML(
					"[R][H]", "Eg2");
		} else {
			int rCount = PolymerisationUtils.countNumberofFragmentPoints(eg2,
					"[R]");
			if (rCount == 1) {
				this.endGroup2 = PolymerisationUtils.convertSmileString2CML(
						eg2, "Eg2");
			} else if (rCount != 1) {
				throw new IllegalArgumentException(
						"End groups are only allowed one R group!");
			}
		}

	}

	/**
	 * Checks CMLMolecule set up is ok for constructor.
	 * 
	 * @param repeatUnitAValue
	 *            - CMLMolecule representation of repeat unit A.
	 * @param repeatUnitBValue
	 *            - CMLMolecule representation of repeat unit B.
	 * @param endGroup1Value
	 *            - CMLMolecule representation of end group 1.
	 * @param endGroup2Value
	 *            - CMLMolecule representation of end group 2.
	 */
	private void checkCMLMoleculeSetUp(final CMLMolecule repeatUnitAValue,
			final CMLMolecule repeatUnitBValue,
			final CMLMolecule endGroup1Value, final CMLMolecule endGroup2Value) {
		if (repeatUnitAValue == null) {
			throw new NullPointerException(
					"No repeat unit entered, please enter a repeat unit");
		} else {
			if (PolymerisationUtils.getRElements(repeatUnitAValue).size() == 2) {
				this.repeatUnitA = PolymerisationUtils.addId(repeatUnitAValue,
						"RUA");
			} else {
				throw new IllegalArgumentException(
						"Repeat unit needs two R groups!");
			}
		}
		if (repeatUnitBValue == null) {
			throw new NullPointerException(
					"No repeat unit entered, please enter a repeat unit");
		} else {
			if (PolymerisationUtils.getRElements(repeatUnitBValue).size() == 2) {
				this.repeatUnitB = PolymerisationUtils.addId(repeatUnitBValue,
						"RUB");
			} else {
				throw new IllegalArgumentException(
						"Repeat unit needs two R groups!");
			}
		}
		if (endGroup1Value == null) {
			this.endGroup1 = PolymerisationUtils.convertSmileString2CML(
					"[R][H]", "Eg1");
		} else {
			if (PolymerisationUtils.getRElements(endGroup1Value).size() == 1) {
				this.endGroup1 = PolymerisationUtils.addId(endGroup1Value,
						"Eg1");
			} else {
				throw new IllegalArgumentException(
						"End groups need to have 1 R group!");
			}
		}
		if (endGroup2Value == null) {
			this.endGroup2 = PolymerisationUtils.convertSmileString2CML(
					"[R][H]", "Eg2");
		} else {
			if (PolymerisationUtils.getRElements(endGroup2Value).size() == 1) {
				this.endGroup2 = PolymerisationUtils.addId(endGroup2Value,
						"Eg2");
			} else {
				throw new IllegalArgumentException(
						"End groups need to have 1 R group!");
			}
		}

	}

	/**
	 * Method returns a fully atomistic CML Document representation of the
	 * polymer.
	 * 
	 * @throws UnbuildableException
	 * @return Document.
	 */
	public final Document createCML() throws UnbuildableException {
		Builder builder = new Builder();
		Document doc = null;
		try {
			doc = builder.build(randomCopolymer.toXML(), "");
		} catch (ValidityException e) {
			e.printStackTrace();
		} catch (ParsingException e) {
			e.printStackTrace();
		} catch (IOException e) {
			e.printStackTrace();
		}
		return doc;
	}

	/**
	 * Creates CMLFragment used to create PML documents.
	 * 
	 * @return CMLFragment
	 */
	public CMLFragment getFragment() {
		ArrayList<CMLAtom> eg1R = PolymerisationUtils.getRElements(endGroup1);
		ArrayList<CMLAtom> eg2R = PolymerisationUtils.getRElements(endGroup2);
		CMLFragment parentFragment = new CMLFragment();
		// Get count of RUA and RUB
		int repeatUnitACount = 0;
		int repeatUnitBCount = 0;
		for (int i = 0; i < ids.size(); i++) {
			if (ids.get(i) == "RuA") {
				repeatUnitACount++;
			}
			if (ids.get(i) == "RuB") {
				repeatUnitBCount++;
			}
		}
		// Create molecule arrays based on count
		ArrayList<CMLMolecule> repeatUnitAMolecules = PolymerisationUtils
				.getRepeatUnitMolecules(repeatUnitACount, repeatUnitA, "RuA");
		ArrayList<CMLMolecule> repeatUnitBMolecules = PolymerisationUtils
				.getRepeatUnitMolecules(repeatUnitBCount, repeatUnitB, "RuB");
		ArrayList<ArrayList<String>> repeatUnitARGroups = PolymerisationUtils
				.getRgroupConnectors(repeatUnitAMolecules);
		ArrayList<ArrayList<String>> repeatUnitBRGroups = PolymerisationUtils
				.getRgroupConnectors(repeatUnitBMolecules);
		// Create blank fragments
		// End group 1
		CMLMolecule endGroup1Blank = new CMLMolecule();
		endGroup1Blank.setRef("Eg1");
		CMLFragment endGroup1BlankFragment = new CMLFragment();
		endGroup1BlankFragment.addMolecule(endGroup1Blank);
		CMLFragment endGroup2BlankFragment = new CMLFragment();
		// Repeat unit A
		ArrayList<CMLMolecule> blankRepeatUnitsA = new ArrayList<CMLMolecule>();
		ArrayList<CMLFragment> blankFragmentsA = new ArrayList<CMLFragment>();
		// Repeat unit B
		ArrayList<CMLMolecule> blankRepeatUnitsB = new ArrayList<CMLMolecule>();
		ArrayList<CMLFragment> blankFragmentsB = new ArrayList<CMLFragment>();
		PolymerisationUtils.addBlankFragments(repeatUnitAMolecules,
				blankRepeatUnitsA, repeatUnitBMolecules, blankRepeatUnitsB,
				blankFragmentsA, blankFragmentsB);
		// End group 2
		CMLMolecule endGroup2Blank = new CMLMolecule();
		endGroup2Blank.setRef("Eg2");
		endGroup2BlankFragment.addMolecule(endGroup2Blank);
		// Join end group 1 to repeat unit
		String moleculeRefs = "PREVIOUS NEXT";
		CMLJoin endGroup1ToRepeatUnit = new CMLJoin();
		String[] eg1JoinValues = new String[2];
		eg1JoinValues[0] = eg1R.get(0).getId();
		eg1JoinValues[1] = repeatUnitARGroups.get(0).get(0);
		endGroup1ToRepeatUnit.setAtomRefs2(eg1JoinValues);
		endGroup1ToRepeatUnit.processMoleculeRefs2AndAtomRefs2(endGroup1Blank,
				repeatUnitAMolecules.get(0));
		endGroup1ToRepeatUnit.setMoleculeRefs2(moleculeRefs);
		parentFragment.appendChild(endGroup1BlankFragment);
		parentFragment.addJoin(endGroup1ToRepeatUnit);
		parentFragment.appendChild(blankFragmentsA.get(0));
		// Join up molecules in sequence [fragment -> join -> fragment]
		int repeatUnitAadded = 0;
		int repeatUnitBadded = -1;
		for (int i = 1; i < ids.size(); i++) {
			if (ids.get(i) == "RuA" && ids.get(i - 1) == "RuA") {
				CMLJoin repeatUnitJoin = new CMLJoin();
				repeatUnitJoin.setAtomRefs2(repeatUnitARGroups.get(
						repeatUnitAadded).get(1)
						+ " "
						+ repeatUnitARGroups.get(repeatUnitAadded + 1).get(0));
				repeatUnitJoin.processMoleculeRefs2AndAtomRefs2(
						repeatUnitAMolecules.get(repeatUnitAadded),
						repeatUnitAMolecules.get(repeatUnitAadded + 1));
				repeatUnitJoin.setMoleculeRefs2(moleculeRefs);
				parentFragment.addJoin(repeatUnitJoin);
				repeatUnitAadded++;
				parentFragment.appendChild(blankFragmentsA
						.get(repeatUnitAadded));
			} else if (ids.get(i) == "RuB" && ids.get(i - 1) == "RuA") {
				CMLJoin repeatUnitJoin = new CMLJoin();
				repeatUnitJoin.setAtomRefs2(repeatUnitARGroups.get(
						repeatUnitAadded).get(1)
						+ " "
						+ repeatUnitBRGroups.get(repeatUnitBadded + 1).get(0));
				repeatUnitJoin.processMoleculeRefs2AndAtomRefs2(
						repeatUnitAMolecules.get(repeatUnitAadded),
						repeatUnitBMolecules.get(repeatUnitBadded + 1));
				repeatUnitJoin.setMoleculeRefs2(moleculeRefs);
				parentFragment.addJoin(repeatUnitJoin);
				repeatUnitBadded++;
				parentFragment.appendChild(blankFragmentsB
						.get(repeatUnitBadded));
			} else if (ids.get(i) == "RuB" && ids.get(i - 1) == "RuB") {
				CMLJoin repeatUnitJoin = new CMLJoin();
				repeatUnitJoin.setAtomRefs2(repeatUnitBRGroups.get(
						repeatUnitBadded).get(1)
						+ " "
						+ repeatUnitBRGroups.get(repeatUnitBadded + 1).get(0));
				repeatUnitJoin.processMoleculeRefs2AndAtomRefs2(
						repeatUnitBMolecules.get(repeatUnitBadded),
						repeatUnitBMolecules.get(repeatUnitBadded + 1));
				repeatUnitJoin.setMoleculeRefs2(moleculeRefs);
				parentFragment.addJoin(repeatUnitJoin);
				repeatUnitBadded++;
				parentFragment.appendChild(blankFragmentsB
						.get(repeatUnitBadded));
			} else if (ids.get(i) == "RuA" && ids.get(i - 1) == "RuB") {
				CMLJoin repeatUnitJoin = new CMLJoin();
				repeatUnitJoin.setAtomRefs2(repeatUnitBRGroups.get(
						repeatUnitBadded).get(1)
						+ " "
						+ repeatUnitARGroups.get(repeatUnitAadded + 1).get(0));
				repeatUnitJoin.processMoleculeRefs2AndAtomRefs2(
						repeatUnitBMolecules.get(repeatUnitBadded),
						repeatUnitAMolecules.get(repeatUnitAadded + 1));
				repeatUnitJoin.setMoleculeRefs2(moleculeRefs);
				parentFragment.addJoin(repeatUnitJoin);
				repeatUnitAadded++;
				parentFragment.appendChild(blankFragmentsA
						.get(repeatUnitAadded));
			} else if (ids.get(i) == "Eg2" && ids.get(i - 1) == "RuA") {
				CMLJoin repeatUnitToEndGroup = new CMLJoin();
				repeatUnitToEndGroup.setAtomRefs2(repeatUnitARGroups.get(
						repeatUnitAadded).get(1)
						+ " " + eg2R.get(0).getId());
				repeatUnitToEndGroup.processMoleculeRefs2AndAtomRefs2(
						repeatUnitAMolecules.get(repeatUnitAadded),
						endGroup2Blank);
				repeatUnitToEndGroup.setMoleculeRefs2(moleculeRefs);
				parentFragment.addJoin(repeatUnitToEndGroup);
				parentFragment.appendChild(endGroup2BlankFragment);
			} else if (ids.get(i) == "Eg2" && ids.get(i - 1) == "RuB") {
				CMLJoin repeatUnitToEndGroup = new CMLJoin();
				repeatUnitToEndGroup.setAtomRefs2(repeatUnitBRGroups.get(
						repeatUnitBadded).get(1)
						+ " " + eg2R.get(0).getId());
				repeatUnitToEndGroup.processMoleculeRefs2AndAtomRefs2(
						repeatUnitBMolecules.get(repeatUnitBadded),
						endGroup2Blank);
				repeatUnitToEndGroup.setMoleculeRefs2(moleculeRefs);
				parentFragment.addJoin(repeatUnitToEndGroup);
				parentFragment.appendChild(endGroup2BlankFragment);
			}
		}
		return parentFragment;
	}

	/**
	 * Creates a list of fragments to build the polymer.
	 */
	private void setUpComponents() {
		ids = randomCopolymerSeparateComponents();
	}

}
